/*******************************************************************************
 * DANIEL'S ALGORITHM IMPLEMENTAIONS
 *
 *  /\  |  _   _  ._ o _|_ |_  ._ _   _ 
 * /--\ | (_| (_) |  |  |_ | | | | | _> 
 *         _|                      
 *
 * INTEGER -- OF ARBITARY LENGTH
 *
 * Retrieved from: 
 * http://en.literateprograms.org/Arbitrary-precision_integer_arithmetic_(C)?oldid=16902
 *
 ******************************************************************************/

#ifndef __INTEGER_H__
#define __INTEGER_H__

#include <stdlib.h>  
#include <stdint.h>
#include <stdbool.h>
#include <string.h>  
#include <limits.h>
#include <math.h>
#include "generic.h"

namespace alg
{
	class Integer
	{
	private:
		typedef unsigned short component_t;
		typedef unsigned long double_component_t;

		#define Max_COMPONENT  ((component_t)(-1))
		#define COMPONENT_BITS  (sizeof(component_t)*CHAR_BIT)

		#define LOG_2_10    3.3219280948873623478703194294894

		component_t* c;    /* least-significant word first */
		int num_components;

	public:
		Integer(int components) 
		{
			num_components = components;
			c = new component_t[components];
			memset(c, 0, sizeof(component_t)*num_components);
		}

		Integer(const Integer & rhs)
		{
			num_components = rhs.size();
			c = new component_t[num_components];
			memcpy(c, rhs.components(), num_components*sizeof(component_t));
		}

		~Integer() 
		{
			delete [] c;
		}

		inline const component_t & operator[] (int i) const { return c[i]; }
		inline component_t & operator[] (int i) { return c[i]; }

		inline const component_t * components() const { return c; } 
		inline uint32_t size() const { return num_components; }
		
		static const Integer from_string(const char* s) 
		{
			Integer result((int)ceil(LOG_2_10*strlen(s)/COMPONENT_BITS));
			Integer digit(1);

			int i;
			for (i = 0; s[i] != '\0'; i++) {
				result = result*10;
				digit[0] = s[i] - '0';
				result = result + digit;
			}

			return result;
		}

		char * to_string() const
		{
			Integer x = (*this);
			int i, result_len;
			char* result = new char[(int)ceil(COMPONENT_BITS*size()/LOG_2_10) + 2];

			Integer ten(1);
			ten[0] = 10;

			if (x.is_zero()) {
				strcpy(result, "0");
			} else {
				for (i = 0; !x.is_zero(); i++) {
					result[i] = (char)(x%10) + '0';
					x=x/10;
				}
				result[i] = '\0';
			}

			result_len = strlen(result);
			for(i=0; i < result_len/2; i++) {
				char temp = result[i];
				result[i] = result[result_len - i - 1];
				result[result_len - i - 1] = temp;
			}

			return result;
		}

		bool is_zero() 
		{
			uint32_t i;
			for(i=0; i < size(); i++) {
				if ((*this)[i] != 0) return false;
			}
			return true;
		}

		// Integer Assignment
		Integer & operator= (const Integer & source)
		{
			memmove(c, source.components(), sizeof(component_t)*Min(source.size(), size()));

			if (size() > source.size()) {
				memset(c + source.size(), 0, sizeof(component_t)*(size() - source.size()));
			}

			return (*this);
		}

		const Integer operator+ (const Integer & rhs)
		{
			Integer result(Max(size(), rhs.size())+1);

			double_component_t carry = 0;
			uint32_t i;
			for(i=0; i<size() || i<rhs.size() || carry != 0; i++) {
				double_component_t partial_sum = carry;
				carry = 0;
				if (i < size())  partial_sum += (*this)[i];
				if (i < rhs.size()) partial_sum += rhs[i];
				if (partial_sum > Max_COMPONENT) {
					partial_sum &= Max_COMPONENT;
					carry = 1;
				}
				result[i] = (component_t)partial_sum;
			}
			for ( ; i < result.size(); i++) { result[i] = 0; }

			return result;
		}

		const Integer operator- (const Integer & right)
		{
			Integer result(Max(size(), right.size()));

			int borrow = 0;
			uint32_t i;
			for(i=0; i<size(); i++) {
				double_component_t lhs = (*this)[i];
				double_component_t rhs = (i < right.size()) ? right[i] : 0;
				if (borrow) {
				if (lhs <= rhs) {
					/* leave borrow set to 1 */
					lhs += (Max_COMPONENT + 1) - 1;
				} else {
					borrow = 0;
					lhs--;
				}
			}
				if (lhs < rhs) {
					borrow = 1;
					lhs += Max_COMPONENT + 1;
				}
				result[i] = lhs - rhs;
			}
			for ( ; i < result.size(); i++) { result[i] = 0; }

			return result;
		}

		const Integer operator* (const component_t & rhs)
		{
			Integer result(size()+1);

			double_component_t carry = 0;
			uint32_t i;
			for(i=0; i<size() || carry != 0; i++) {
				double_component_t partial_sum = carry;
				carry = 0;
				if (i < size())  partial_sum += (*this)[i]*rhs;
				carry = partial_sum >> COMPONENT_BITS;
				result[i] = (component_t)(partial_sum & Max_COMPONENT);
			}
			return result;
		}

		const Integer operator* (const Integer & rhs)
		{
			Integer result(Max(size(), rhs.size())*2);

			uint32_t i, lidx, ridx;
			double_component_t carry = 0;
			uint32_t max_size_no_carry;
			uint32_t left_max_component  = size() - 1;
			uint32_t right_max_component = rhs.size() - 1;
			while((*this)[left_max_component] == 0) left_max_component--;
			while(rhs[right_max_component] == 0) right_max_component--;
			max_size_no_carry = left_max_component + right_max_component;
			for(i=0; i <= max_size_no_carry || carry != 0; i++) {
				double_component_t partial_sum = carry;
				carry = 0;
				lidx = Min(i, left_max_component);
				ridx = i - lidx;
				while(lidx >= 0 && ridx <= right_max_component) {
					partial_sum += ((double_component_t)(*this)[lidx])*rhs[ridx];
					carry += partial_sum >> COMPONENT_BITS;
					partial_sum &= Max_COMPONENT;
					lidx--; ridx++;
				}
				result[i] = partial_sum;
			}
			for ( ; i < result.size(); i++) { result[i] = 0; }
			return result;
		}

		const Integer operator/ (component_t rhs) 
		{
			Integer result(size());
			double_component_t dividend = 0;
			int i;
			for (i = size() - 1; i >= 0; i--) {
				dividend |= (*this)[i];
				result[i] = dividend/rhs;
				dividend = (dividend % rhs) << COMPONENT_BITS;
			}

			return result;
		}

		component_t operator% (component_t right) 
		{
			double_component_t mod_two_power = 1;
			double_component_t result = 0;
			uint32_t i, bit;
			for(i=0; i<size(); i++) {
				for(bit=0; bit<COMPONENT_BITS; bit++) {
					if (((*this)[i] & (1 << bit)) != 0) {
						result += mod_two_power;
						if (result >= right) {
							result -= right;
						}
					}
					mod_two_power <<= 1;
					if (mod_two_power >= right) {
						mod_two_power -= right;
					}
				}
			}
			return (component_t)result;
		}

		const Integer operator% (const Integer & rhs) 
		{
			Integer result = rhs;
			Integer mod_two_power(rhs.size() + 1);

			uint32_t i, bit;
			mod_two_power[0] = 1;
			for(i=0; i<size(); i++) {
				for(bit=0; bit<COMPONENT_BITS; bit++) {
					if (((*this)[i] & (1 << bit)) != 0) {
						result = mod_two_power+result;
						if (result.compare(rhs) >= 0) {
							result = result - rhs;
						}
					}
					mod_two_power.shift_left_one_integer();
					if (mod_two_power.compare(rhs) >= 0) {
						mod_two_power = mod_two_power - rhs;
					}
				}
			}

			return result;
		}

		int compare(const Integer & rhs)
		{
			uint32_t i = Max(size() - 1, rhs.size() - 1);
			for ( ; i >= 0; i--) {
				component_t left_comp =
					(i < size()) ? (*this)[i] : 0;
				component_t right_comp =
					(i < rhs.size()) ? rhs[i] : 0;
				if (left_comp < right_comp)
					return -1;
				else if (left_comp > right_comp)
					return  1;
			}
			return 0;
		}

private:
		void shift_left_one_integer() 
		{
			uint32_t i;
			(*this)[size() - 1] <<= 1;
			for (i = size() - 2; i >= 0; i--) {
				(*this)[i + 1] |= (*this)[i] >> (COMPONENT_BITS - 1);
				(*this)[i] <<= 1;
			}
		}

		void shift_right_one_integer() 
		{
			uint32_t i;
			(*this)[0] >>= 1;
			for (i = 1; i < size(); i++) {
				(*this)[i - 1] |= ((*this)[i] & 1) << (COMPONENT_BITS - 1);
				(*this)[i] >>= 1;
			}
		}
	};
}

#endif //
